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ONE-DIMENSIONAL TWO-PHASE FLOW IN COMBUSTION CHAMBER 
OF SOLID PROPELLANT ROCKET MOTORS 

CHANG XIANQI 

ABSTRACT 

'■In this paper, a numerical solution of basic equation for 
one dimensional two-phase nonequilibrium flow in a combustion 
chamber of solid propellant rocket motors is discussed in detail, 

combust ion 

the effect of particle size on flow field in chamber and pressure 

time curves is analyzed, and some useful conclusions are obtained 

in comparison with results of one dimensional two-phase constant 
combustion 

lag flow in.chamber. It is useful for predicting pressure-time 
/» 

curves accurately and providing accurate boundary conditions 
for the calculation of two-phase flow through the nozzle. ( ii . • 











Nomenclature 


/II 


A - duct cross section area 
- Burning area of charge 
Aq - Thermal equivalent of work 
A i - Area of the nozzle throat 
M - Mean molecular weight of gas 
n - Pressure exponent 
P - Pressure 
P s - Total pressure 
b * Burning speed coefficient 
C* - Characteristic speed of propellant 

Cj - Particulate specific heat 

Cp g - Specific heat of gaseous phase at constant pressure 
g - Gravitational acceleration 
h - Enthalpy per unit mass 

H s - Total enthalpy in 1 kg of two-phase mixture 
H sg - Total enthalpy in 1 kg of gas 
H sv - Total enthalpy in 1 kg of liquid 
k - Ratio of specific heat of gas 
K - Particle velocity lag coefficient 

vk 

K i - Surface - throat ratio J\ = - 

At 

1 - Charge length 

L- Particle temperature lag coefficient 
m - Mass flow rate of two-phase mixture 
rrig - Mass flow rate of gas 
q - Heat flux of particles per unit mass 

r - Burning rate 
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r p - Radius of particles 

R g - Gas constant of gaseous phase 

S - Duct circumference 

t - Time 

T - Temperature 

To ~ Total temperature of charge head gaseous phase 
T s - Total temperature of gaseous phase 
v- Velocity 
x - Axial coordinate 
X - Particle resistance per unit mass 
y c - Outer radius of charge 
p - Density 

p t - Density of propellant 
m. p - Mass flow rate of condensed phase 
M g - Mach number of gaseous phase 
Pmp ‘ Density of AL 2 0 3 material 

£ - Fraction of particle mass flow rate, e = 

m 

Ay - Coefficient of thermal conductivity of gaseous phase 
My ‘ Coefficient of dynamic viscosity of gaseous phase 

SUBSCRIPTS 

- Gas phase 

o- Cross section of charge head 
t- Cross section of nozzle throat 
P - Condensed piiase 

L - Cross section at the charge tip 

i - initial value 
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1. PREFACE 


Modern composite solid propellants are mixed with a certain 
amount of aluminum powder for increasing energy and decreasing 
unsteady burning. When the aluminum mixed propellant burns, 

Al 2 0 3 particles are formed in liquid phase. The weight fraction 
can be up to 30-40%. Therefore, in combustion chambers and jet 
nozzles, the products of combustion actually are a mixture of 
gas and liquid. 

There have been a lot of publications regarding two-phase 
flow in jet nozzles of a rocket motor. Most of them are devoted 
to one dimensional two-phase flow and two dimensional axial sym¬ 
metrical two-phase flow in nozzles. There has been some progress. 
Two-phase flow in a combustion chamber is involved with mass, 
and has some new characteristics. The author has done research 
on one dimensional two-phase constant lag flow in a combustion 
chamber [4] , analyzing the particle velocity lag effects on the 
combustion chamber processes. This article studies the one 
dimensional two-phase nonaquilibrium flow in a combustion chamber. 
Based upon the principle equations of one dimensional mass involving 
two-phase flow in a combustion chamber., it discusses the numerical 
solutions of the equations in detail, analyzes particle size 
effect on the chamber flow field and the pressure-time curve, 
and compares the results with that of constant lag flow. 

Some practical conclusions are reached. This gives precision 
in predicting pressure-time curve, and more accurate boundary 
conditions for calculation of two-phase flow in jet nozzles. 
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II. FUNDAMENTAL EQUATIONS 


Assume: 
1. Flow 


3. 


is one-dimensional and steady-state; 

2. The friction and heat loss on the duct wall are negligible; 

A^O^ particles are spherical, uniform, and in the liquid 
phase; the particle volume and the Brownian motion effects 
on pressure are negligible; 

The gas phase is an ideal gas at freezing point except in 
contact with particles and otherwise is inviscid; 

No mass exchange between two phases; 

The specific heat of gas and particles i s constant /13 


Based upon the above assumptions, the fundamental equations 
of two-phase flow in the combustion chamber are obtained 
as follows DJ and C3j: 

Gas phase - 
Mass Equation: 

d . 

( 1 ) 


— e)p T rS 


Momentum equation: 


£(P*vlA)--A^-Xp,A 


( 2 ) 


Energy equation: 

d / u* \ 1 

p t gVfA[hg + Aq { = (1 —e)p T rSgHsg — AoXppVpA + qprA 


dx\ 

in which 


(3) 


Condensed phase- 
Mass Equation: 


r = bP" 

ho — C n T t 


-j^iprVpA)** ep T rS 


(4) 






\ s»» 
















Momentum equation: 




Energy equation: 

Af prgvrA(hr+A<) eprrSgHs? + AqX pt-vpA — qppA 


dx ^ rsf -r.. V r . -V 2g 

Two-phase mixture- 
Mass equation: 


dth c 

jj-e.rS 


Momentum equation: 


dP 

-j-iPiv'.A +p,v',A)~ -A-fe 


Energy equation: 


&*• (‘- + ^ fj) £)]-' , ' rSHj 


in which 


»*»* = P&VfA 
thp = PrVpA 
fh ” 1 fftf + 


III. SOME SUPPLEMENTARY RELATIONSHIPS 

Equations (1) — (6) &re not closed, therefore, 
we introduce the following supplementary relationships: 


1 . 


Gas-phase condition equation 

for an ideal gas — P= gR>T 
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2 . 


Condensed-phase condition equation. 

When the Al^O^ particle temperature is greater than the melting 
point (T Pm — 23l8*K) assuming its specific heat is constant, 
then 

hr m *h Pm + Ci(T P — T Pm ) (11) 

in which 

- Specific heat of liquid A^O^ particle, 0.34327 KCal/kg. 


degree K; h - Enthalpy of liquid A^O^ at T pm , 


pm' 


876.9498 


KCal/kg 


From (11) we have 

dh P = C t dT P (12) 

The particle resistance X[l,2] per unit mass 

X=‘A p (v, — v t ) (13) 

Under the condition of combustion chambers, the particles 
carryout Stokes flow, at this moment 


A r = - -&*- 


2 r' P p mP 

2 

in which p, = l .208 x Kr , 7VA7 , - , (kg, sec/m ) 

4. The heat flux q [1,2] of particles per unit mass 

q=B P (T'-T.) (14) 

in which B P =^—— 

5. The total enthalpy H s in 1 kq of two-phase mixture H s = (\ — e)H Sg + eH SP 
From the assumption 2, H g is constant along the channel, 

therefore, it is convenient to use the parameters 0» = i>= 0. T e = T P = T % ) 
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V 

V 

* 

X 


at the 


charge 


origin (x=0), to express H . Therefore 


H s — (\ — e)C Pf T, + e[h Pm + C i{T, T Pm )] 


(15) 


a 


3 


V 

3 
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IV. COMPUTATIONAL EQUATIONS 


Let 


/Cs-' (OsSATCl) 

(16) 

v t 


(o<l<d 

* •“ 1 9 

(17) 


Therefore particle velocity lag ^ = 1 —/C 

particle temperature -lag = 

T,—T t 

Here k,L are defined as particle velocity lag coefficient 
and temperature lag coefficient respectively. 

Ignoring the effect of burning erosion, assuming the cross 
section area A of the charge duct is constant along 
the longitudinal axis, after an elaborate manipulation, the 
following numerical solution is obtained based upon the 







in which 


H s , = C„Ts 


The system of equation (18) contains 3 ordinary differential 
equations and 8 algebraic equations, with the unknowns.^ ^ r ' v *' 


p»*P+ Ps , T s% MgRih, un( j er t h e given boundary conditions, which 


can be solved by the Runge-Kutta numerical method. 


V. INITIAL CONDITION AND BOUNDARY CONDITIONS 


The initial condition of the equation: 


When t=0, S=S., A=A. . The value of A., S- can be obtained by 
i' 10 11 1 

the actual shape of a charge duct, for a circular 


cross section of the duct. 



Si— 2xy t 
Ai—vy', 


y^ is the duct initial radius. 


The boundary conditions of equations: 


If the charge 


fills the annular space between two closed 


end coaxial cylinders, the boundary conditions are (refer to 


Figure 1); at the charge head (x=0) cross section 


tV# = i’*, = 0 


fh, = 0 


Tg t =T Ft =T 


P = P. 


gR,T, 


p„ = _£_ Pm 

1 ~e v r , 1 — eK, 






v_‘0 








I 

I 

I 

in which, the total temperature of the charge head gaseous phase Tq is equal 
to the burning temperature of the charge.. P q is the correspond¬ 

ing pressure. Before solving the equation, it is unknown. Therefore, 
its value has to be determined by the iteration process in the 
numerical solution. The particle density Pr * at the charge 
head is ^ an uncertain value. To determine this value, the 
value of K q has to be obtained beforehand. 



At the charge end (x*l) cross section the mass flow rate m of the 
two-phase mixture passing through the cross section area at the 
charge end (x-1), should be equal to the mass 

flow rate m of the two-phase mixture passing through the throat 
of the jet nozzle, that is 

ih L =‘tti, ( 20 ) 

ft» L can be determined by the parameters of 
of the x=I cross section, that is 

( 21 ) 

Ttl l = A( pg lV g L + P r L*-' p L ) 
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The present article is mainly about two-phase nonequilibrium 


flow in a combustion chamber. To avoid the numerical solution 
for one dimensional nonequilibrium two-phase flow in the jet 
nozzle, in the determination of the boundary condition 
of Che charge end, it is assumed that the 

flow in the jet nozzle is a one dimensional two-phase constant 
lag flow, therefore: 

gR,T,C {22) 


in which 


C ~ 1 + iZr e {K C*( 1 - K ) + AT] + (k- 1) 6LD} 


D = 


1 + 


1 —e 


1+-J-SL 
1 — e 



f _l + (fe-l)g 





where the value of K, L should be taken as the corresponding 
value at the end of the field length. 


VI. THE DETERMINATION OF K , L 
( ° ° 

For obtaining the distribution of particle speed lag and 
temperature lag along the duct,, and the numerical solution 

for a system of equations (18), the values of K, L at the 









charge head (x=0) have to be determined. 

However, at x-0, both K and L are indeterminate forms, therefore, 
it is necessary to consider that at the charge head 

i 

) 

l 

! 

I 

1. The determination of K 

o 

' from 1 1 Hopital's rule 


K» 




(a) 


After transforming the first equation in the system of equations 
(18) , we find that , at x-*0* \ is also 

an indeterminate form of So 


Vr\ m * t +-Ar —-;- *>f 


For convenience, the subscript *-»0 + has been omitted in 
the following expression. After rearrangement it reads 


I 



- v;=*o 


(b) 


With the help of the boundary condition (19), from the second 
equation in the system of equations (18) the following is 
derived that at 


in which r.^bPl 


P..A 


(c) 


Substitute (c) into (b), it is obtained that at 


x-»0 + 






Because in a combustion chamber the flow involves 


an additional mass of particles, if the speed of a particle 
increases, then therefore, (d) should be positive 

which results in the sign in front of the square root being 


Substitute (d), 


then obtain 


(c) into (a) , and set 

^ 4 (l-t)prr,5 

a Ar P»*A 


K.~ 


-l + Jl + 2a 
a 


c 


c 


The determination of L : 

o 

Same as above, from 1' Hopital rule 



It can be derived from the fundamental equations (1)—(6) 


4Tjl 

dx 


T' 1 In —h 4 -A ^ \ r — Tt) 

. T ^^ Hsr -h, + A Q ~) VpgCi 


herein, T'A*—* is also an indeterminate form of 


0 

0. 


From 1'Hopital rule it is obtained as 


T‘ r \~,+ =-T’r 


B r T f 1»—» * T , \ fi 4 
gCi »*!.»•♦ 


Omitting the subscript , then 













from equation (f) it is obtained 


I 



(25) 


With the expressions of equations (d) and (23), the following 


is obtained: 


. * B, 1 

P gC, A , ( — 1 + %/ 1 + 2a) 


(26) 


VII. NUMERICAL SOLUTION OF THE SYSTEM OF EQUATIONS 

Under the boundary conditions (19) , (20) and the initial 
condition, but use of the fourth order Runge-Kutta method the 
solution for the system of equations (18) can be obtained as 
follows: 

1. At a fixed time t, the gas flow parameters (t>», T,, p,, P, 

Pi, T s , M,),. , the particle flow parameters p and 

the particle lag factor K, L distribution along the 
duct. 

> 

f ... 

f 2. The pressure variation with respect to time in the combustion 

I 

' chamber. 

In the boundary condition (19), the head pressure is unknown apriori, 
and has to be determined with the numerical solution process. The 
procedure is the following: 







(1) Form pi 1 *— (C*p r bK to compute the first approximate value; 

(2) With the value of P ^ solve for the system of equations 

(18), and obtain the variation of parameters of the gas and the 
particle along the x-axis, then form equation (21), (22) and 

solve for m^ and m fc ; 

(3) Let Artiom, — rhi and make a judgment on 

|45|<,, (27) 

I Ihi I 


( is a given allowable error, for instance £^=0.01) to see 
whether it is satisfied. * 

(4) If equation (27) is not true, and a*»»> 0! , take AP 

(AP is the given increase in pressure, for example 

A a p =1 kg/cm x ) , otherwise take p‘»>*p<'>+ AP , repeat the pro¬ 

cedure for solving the system of equations (18) , if the mass 
flow rate satisfies (27), then P ^ is the true value of the 
pressure at the charge head P Q ; 

(5) If the condition (27) still can not be satisfied, the following 
interpolation equation can be used for computing P ^ : 


p< 

* v 


< •>. 


. .> 


then solve for the system of equations (18) , till the mass flow 
rate satisfies equation (27) . Figure 2 shows the flow chart 


of the numerical solution process, in which ‘T* 


cp 


burning speed along the length of the duct. 


is the average 








lV- 







Figure 2. 

Key: (1) start;(2) input;(3) compute P Q ; (4) furnish 
boundary conditions; (5) Runge-Kutta method for solving 
the system of equations (18), (6) output, (7) adjust 
the head pressure (8) No, (9) Yes, 

(10) Yes, (11) No, (12) Stop, (13) Flow Chart of 
Computation. 
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VIII. THE EFFECT OF THE PARTICLE SIZE ON 
INTERNAL BALLISTIC PROPERTIES 

As to a certain solid propellant rocket motor ( e =0.26), 
the pressure-time curve and the flow field in the combustion 
chamber have been computed for different sizes of particles. 

Figure 3 shows the effect of the particle size on the Pressure- 
Time curve. It can be observed in the figure that the pressure 


in the combustion chamber decreases as the particle size increases. 

At t=0 sec., the variation of the combustion chamber pressure 

P Q with respect to the particle radius is shown in Figure 4. 

The decrease of the combustion chamber pressure will lead to 

a slowdown of the burning speed, increase the charge burning 

period and decrease the mass flow rate. 

At t=0 sec., the variation of the gas parameters along the 

x-axis for the different particle sizes is shown in Figure 5. 

Mg solely increases with the distance, the rest of the parameters 

decrease with the distance. Because P Q decreases when the particle 

size increases, therefore P,P , p are affected by the particle 

s g 

size distinctively. The larger the particle size, the smaller 
are the values of them. But the particle size effect on Tg is 
negligible. 

Figure 6 shows at t=0 sec. the distribution of the gas speed 
, g and the particle speed Vp along the duct under 
different particle size conditions. Because the flow 
in the combustion chamber is involved with an increase in mass, 
both v-g and increase along the length of the duct. But the particle 


17 





v v 







size effect on and \f is different when particle size increases 

<3 P 

\Tg increases, but decreases. This is because the particle 
of unit mass has a resistance to the gas of 92 p- . When the 

' P. 

particle size increases, X decreases, so the gas speed increases. 

At this moment the torque exerted on the particle by gas decreases. 
This results in the particle acceleration decreasing, therefore 
'C'p becomes small. The density and temperature distribution of 
the particle are shown in Figure 7. It can be observed that 
fip hardly decreases along the duct and T^ basically stays 

constant. 

The effect of the particle size on the value of is 
significant: when particle size increases, increases too. 

But the effect on the value of T^ is negligible. 

At t^=0 sec., the particle speed lag factor K and temperature 
lag factor L distributions along the duct are shown 
in Figure 8. It can be observed that the value of K slightly 
decreases along the distance, but the value of L fundamentally 
stays constant. This indicates that in the combustion chamber, 


the particle speed lag (1— K ) increases slightly along the length of 
the duct, but the temperature lag (1-L) barsitsa’Ily stays constant. 

The smaller is the particle size, the smaller is the variation 
of K value along the distance. (Refer to Table 1.) Therefore 
when the particle size is small, the two-phase flow in the 
combustion chamber can be treated as the constant lag flow motion. 
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Table 1. 



Key: ( l) The radius of the particle r 

It 

(micrometer), ( » The decrease of K at the charge end 

Table 2 lists, for different particle sizes, the comparison between 
the numerical solutions and constant lag flow computations of 
some important parameters in two-phase flow in a combustion 
chamber. Compared with the numerical solution, the combustion 
chamber pressure P q and the total pressure at the end 
charge end obtained from the constant lag flow 

computation are comparatively low, and the rest of the parameters 
are comparatively high, but the deviation is less than 1%. 


Table 2. 



2 2 

Key: (2) Micrometer, (3) Kg/cm , (4) kg/cm , 

(5) m/sec, (6) m/sec, (7) ( 9),( 11),(13) numerical solution, 

(8) (10),(12) (14) constant lag flow. 
















From equation (24), (25) it becomes known the values of K Q , L Q 

are functions of the particle radius r . Their variation with 

respect tor is shown in Figure 9. It can be seen in the Figure 
P 

K , L values decrease when the particle size increases, which 
o o 

means, the particle speed lag and temperature lag increases as 
the particle size increases. 

Besides, K q < L q , in the nonequilibrium flow in the combustion 
chamber, the lag of the particle speed is greater than the lag 
of the particle temperature. 

IX. CONCLUSION 

Through the above discussion the following conclusions are 
therefore obtained: 

1. Because of the effect of the two-phase flow, the pressure 

of the combustion chamber decreases, and the larger the particle 
size, the smaller is the combustion chamber pressure. 

2. The two-phase flow affects the flow field in the combustion 
chamber greatly. When the particle size increases, the gas 
speed increases, the pressure, total pressure and the gas 
density all decrease. When the particle speed decreases, 
the density increases. The effect of the particle size on 
the temperature of both gas and particle is very 

small. 

3. In the combustion chamber the lag of the particle speed increases 
along the length of the duct, the lag of the temperature basically 
keeps constant. Besides, the particle speed lag is greater 
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SUMMARY 


On the basis of the evidence presented, a steady state combustion 
model (l) of solid state propellents with negative pressure coefficient 
burn speed characteristics allows us to derive a new pressure response 
function formula. This can be used to explain pressure pairing phenomena 
for propellents whose combustion speed pressure exponents are zero, 
positive, and negative. The burning propellent is divided into two 
sections: the first section is a structure composed of oxidizing, agent 

covered by nolten binding agent and the corresponding binding agent 
surface. The other section is formed by uncovered oxidizer surface and 
remaining binding agent surface. This model is different from the 
various types of models in the past. In the combustion on the surface of 
the first type of section described above, consideration has been given 
to the oxidizer, under conditions in which it is covered by molten 
binding agent, so that it is considered to be in a state of opposed 
gasification and congealed phase reaction. Therefore, the real section 
of the pressure response function which is obtained, when the pressure 
exponent of the propellent steady state fuel speed is zero or has a 
negative value, is also capable of being a positive value. When we made 
use of the expression obtained for the pressure response function in 
experimentation with the propellent (S04-5A)and made qualitative 
calculations, the results of these calculations satisfactorily explain 
the phenomenon of propellents with negative pressure exponents still being 
unstable in coabustion when rust of the surface area of the oxidizer is 
covered by molten binder agent. This not only overcomes the weaknesses 
which all expressions for pressure response functions had in the past 
when used with negative pressure exponent propellents, but also, in a 
way, reflects the accuracy of the solid fuel propellent steady state 
combustion moiel (l) in its combustion speed characteristics for 
positive and negative pressure exponents. 

Explanation of Symbols 


A 


4 
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t . -a • . 


w-V-V. 




















Oxidizer Gasification Rate Indicator Prefactor 


A 


•* 


B 


_ 2R*T l , 
={T.-T,)E/ 


E, Indicator Prefactor in Oxidation Agent Equilibrium Evaporation 
i=>. Te Formula 


w * ) i 
Pexp (*%) 


c opecific Heat 




C + q/E tm 


E,. Energy of Activation for oxidizer Agent Gasification Reaction 

Ef Energy of Activation for the Propellent Gas Phase Combustion 

one s s 

p mQ+ A t (AB — A) 


G Hass fraction of Oxidizer Agent for Completion of Condensation 

t ■ - '.eaction 

K ".peal Constant for Gas Phase Reaction 










m Mass Plow Rate 

„ Combustion Speed Pressure Exponent 
p Pressure 

Qf Unit Mass "Binder Decomposition Reaction Heat 
Q. Boundary Layer Reaction Heat of Unit Mass Oxidizer 
9 Heat of Evaporation for Unit Mass Oxidizer 

<h Heat Plow for Gas Phase Boundary Surface Plow 
Qo Vacant of Heat Released by Unit Mass Gas Phase Reaction 

q< Rate of Release for the Heat of Purification from Boundary Surface 

Reaction 

0. Real Section 

/?' Universal Gas Constant 


Rr Pressure Response Punction 


mlP 


S Ejtuiisji] ( i is the imaginary number unit J —1 ) 

ffl*c 

j Temperature 

T/ Absolute Temperature of the Propellent Flame 
T, Initial Propellent Temperature 


T, D ~opellent Surface Temperature 










W A , Kt) rram molecular weight of oxidizer evaporation gas 

W G Gram Molecular Weight of Gaseous Products from Congealed Phase 
Oxidizer Reaction 


x Distance 

a Mass fraction of Oxidizer in the Propellent 

y Surface Area fraction Composed of Oxidizer Agent Surface Covered by 
Molten Binder Agent 

X Congealed and Solid Phase Meat Conduction Coefficient 
i Mondimensional Distance 


p Congealed and Solid Density 


O Mondimensional Frequency 


co Angular Frequency 


SUPERSCRIPTS 


— Steady State Value or Average Value 


, r, e”tu”baf i)'i Value 


~ Co ml •> < \ anl : tude of Perturbation 




. ' S '. .% /. . 









1 Parameters Corresponding to Area I 


1 Parameters Cor responding to Area II 


3U33CRIPTS 
« imaginary Section 

r Peal Section 

I Parameters Correspond!ng to Area I 

1 Parameters Corresponding to Area II 

s Soundary Surface, Value (Congealed Phase Side) for the Place Where 

x — 0 

1, Introduction 

Sonic instability is the primary form of combustion instability of 
soli I propellent rocket motors. Moreover, this type of sonic instability 
has as its primary source the combustion response of solid fuel 
propellents. Because of this, in the test construction of solid rocket 
actors, the sonic instability which they show is extremely important, and 
the theory and measurements associated with combustion response are I 

indispensable. In linear sonic analysis, on sees shown whether or not a j 

soli 1 -ooket motor will give a timely, consistent pressure perturbation j 

e ol - -ient, forming a sonically unstable fuel. In this type of 

; l tly ' i •, t'l i pressure response function must already be known. 1 

“k 1 o ’ )hl j i ■> * soli i propellent pressure response has received 
I a r.Dint s of research work in both the U.3. and the Soviet Union. 

! 











F.E.C. Culick (2) and N.S. Cohen (3) have already done excellent 
critiques of this problem. Concerning "flatbed" propellent pressure 
response functions, P.E.C. Culick and others (2) (4) have also done large 
amounts of research. However, the author recognizes that the steady 
state combustion model and several assumptions which form the basis of 
the research are still short on experimental data to support them. 
Moreover, several points of theory are still doubtful. Por example, in 
order to explain the phenomenon of sonic vibrations in the combustion of 
flatbed propellents, we take the mass gasification rate for solid phase 
reaction areas and change it for use in the expression below: 

f. 

m = A,P‘>e 

Prom this, we can derive the pressure response function formula: 

r __ nAB + n,(A— 1 ) 

A + ^-(l +A) + AB 


Moreover, when we use this relationship to compare with the results of 
■i <per! ns a ts, we select n, = i.o . However, we know that, when the 

pr i i jure exponent of the lass gasification rate for solid phase reaction 
areas is a positive value, the pressure exponent for its steady state 
combustion speed can hardly be zero. This article does not intend to 
make more criticism of these publications. We will present on the basis 
of a steady state model (1) of solid fuel propellents in terms of 
positive and negative pressure exponents for conbustion speed 
characterist ics, and derive a new pressure response function formula and 
fo-’oe diagram- We will use the i to make a complete explanation of the 
pressure pairing phenoaenon wVc’n exists in various types of solid 
propellents when their combustion speed pressure exponents are zero, 
positive, and negative. 









2, Physical-Chemical i-iodel 

In the same way as was the case with the corresponding steady state 
model (l), we select for use a one dimensional aodel and divide the 
combustion surface of the propellent into two sections. One section is 
exposed of oxidizer agent surface covered by molten binder agent and the 
o v'-'e-sp filing binder agent surface. Let us also stipulate that the 
fraction of surface area which it occupies out of the whole combustion 
surface is v . The other section is then composed of oxidizer surface 
which is not covered and the remaining binder agent surface. These two 
sections are respectively called, for short, Area I and Area II, as shown 
In 7 lg. 1 . At the same time, we take the combustion process for each 
solid propellent area and simplify it into three stages, which 
respectively occur in three different phases (Pig. 2). 

1) The solid phase area on the inside of the solid propellent, to 
which heat is applied, is gasified by the gas phase flame area with the 
congealed phase reaction area supplying the heat. 

2) Within the congealed phase layer between the solid phase and the 
gas phase one sees the development of an oxidizer boundary layer reaction 
and binder agent heat of decomposition which are contained in the 
congealed phase reaction and the gasification reaction. In Area I, 
unlike Area II, one must consider the molten binder agent liquid layer 
which covers the surface of the oxidizer as well as the oxidizer 
underneath It as they exist in the congealed phase reaction related to 
gasification and its opposite process—opposed gasification. 

3) In the gas phase one sees the occurrence of combustion processes 
which include dispersion, mixing, and chemical reaction. 
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Pig. 1 One Dimensional Steady State Model 
1. Dispensed Plame 2. Perturbed Flame 3* Oxydizer Agent 4. Liquid 
Phase 5- Solid Phase 6. Paired Agent Mixed Layer 7. End Plame 8. 

AP Plame 9. Initial Plame 10. Combustion Surface 11. Propellant 12 
Liquid Phase 13* Solid Phase 
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Pig. 2 4 Simplified Two Dimensional Model 1. Congealed Phase Layer 

Solid Phase 4rea of Increased Heat 3* Sas Phase Reaction 4rea 











3, Mathematical Treatment 


For the sake of convenience in the mathematical treatment, we first 
make the following assumptions: 

1 ) The influence of heat damage is negligible. 

2) The solid propellent is noncompressible, of a uniform qualtiy, 
and of the same nature in all directions. 

3) The congealed phase layer between the solid phase and the gas 
phase is an infinitely thin plane. In this plane one sees ocouring the 
oxidizer congealed phase reaction and the gasification reaction (Area I 
•and Area II are different and each has its own rules of reaction) . At 
the same time, one sees occuring the high temperature decomposition of 
binder agent. This plane is called the boundary surface. 

(4) In the three solid, liquid, and gas phases discussed above, the 
specific heat is always a constant. 

5) las phase reactions can be seen to be steady stats processes. 

6) The influence of congealed phase reactions on pressure 
disturbances is negligible. 

7) When the amount of perturbation is sufficiently small, its second 
and higher order small amounts can all be neglected. 

At the same time, in the same fashion as normal models, we choose 
for use a moving coordinate system. These coordinates take the 
instantaneous burn speeds of the propellents and move them into the solid 
phase area. In this way, the origin of the coordinate system can fall 
entirely on the solid-gas phase boundary surface. 

Below we present the respective treatments for Area I and Area II. 

Area I: 

(l) Basic Equations 

1) The unstable thermal conductivity equation in the solid phase 
heat addition area is: 


. fi'-T fiT fiT 

K -W~ mc aT =pc * 


'J ) 
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2) When we consider the existence of molten boundary agent on the 
cohering of oxidizer surface and the oxidizer under the covering as they 
affect congealed phase reactions and gasification reactions, from 
reference article (l), the mass gasification rate of propellent is: 



3) If ve assume a steady state, the mass combustion rate for the gas 
phase reaction area can be expressed as (5): 


,-KP’ex p( 2R*T /) 


4) Energy equilibrium equation for the inside of the infinitely thin 
congealed phase layer: 

The difference between the heat flow <7. from the gas phase 

,dT I 

combustion area and the heat flow A !>x\ transfered into tne 

solid phase should equal the rate of heat release q r , which is produced 
by the reaction in th? congealed phase layer of the congealed phase and 
the gasification reaction, that is, 




Moreover, ;-5 i ag motion coo rd i nates for explanation, Qi = mLQi-c(Tf—T,)2, 
,also, Qc _ m [ c Q I -(i-a)Oflo . In these equations, Q c is the 

amount of heat relent d by the unit a ass gas phase reaction. Q* is 
the algebraic sura of the effective heats of the oxidizer agent congealed 
phase reaction and the gasification reaction. CV is the effective heat 
from the high temperature thermal decomposition of binder igi-it. I ? ■: i 














substitute the equations above, we get 


k-S- +m[aQ.-{\-a)Q F -Q,+c{T l -T,)-\ =0 

OX 


) Linearization Treatment: 

1) ^irst, we carry out a linearization treatment of equation (1 ). 
We postulate: 

m(/) = #ft 4- : 


7\x, 0 = T(x) + /?.{f(>V f *"} ‘ 


If we substitute in equation (l), and eliminate the high order small 
quantities, as well as making use of the steady state condition 


<d'T „dT_ Q 

A r~r — fflc - u 

dx 1 dx 


n, e pint ion .1) changes to become: 




•i i-;e r.'n fa- as - t r e set up for all times ‘ , therefore, one nas 












( 7 ) 


dx' dx “* 


t m!f. 


S-tofe-* 


Then equations (o) and (7) can respectively be written 


dV dt 

d 2 f _df _cf Ma >» dT_ 
dV dt rn di 


2) In the same way, in the case of equation (2), one can set 


m = fn + 7n / = ffl + R,{1he <m< } j 

T.-r. + T\-r. + R.{T.e ,m ') 

P = P + P' = P + R'{Pe‘ m '} ] 


Moreover, one can use the steady state condition 

! ! , / £..J 

m= a T^^' eXp l 1 n /. , g 

V l — G W c ) 

Division with equation (2) gives one 


( 10 ) 


(11 ) 
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1 + — -exp y R , Ti • Tt J 


In the equation 


/ G W 

Cs ^xlIEl w c ) _, 

Peycp(-^r) 


E T' t' 

Also, bee vise of the fact that p i tj r • jr and ' f~ are both 

very small (tP?~ ari< ^ ~R*¥~ are both the or ber of magnitude 
of 10). Moreover, ~~ is then in the range KT'-lO"' • Therefore 
one has 


eX T)f E, . . Tj_ \ ~ | + ^9*, . Tj_ 

exp \ R'T, T.) R'T, T, 

cxp (^?- . i__a__ . 71 

exp \ R'T, T.) R'T, T. 


Substituting in eolation (12) and eliminating the smaller quantities of 
higher order, one, then, obtains 

, „ R’ _ g_ r 

T. h. ? R'T. T. 

m R'T, ' T, C 


As ibove, because of the fact that all of the equivalent forms for time.- 

'ho ill b 1 •• et io, one, lonse q.iently has 


f r 


C+T - T. R'T, Pi 








3) After one carries oat the same type of linearization treatment 
on equation (3), one obtains: 


A P 
t_l a~"P 

r, TT 

2 R*T, 


(H) 


4) In the same way, equation (4) can be changed to become; 




(15) 


(3) Solution of Equations 

Because equation (9) and its boundary condition are: 


d'T dT 
dV di = ° 


t= —OO t ime T = T . 
'£*=0 time 2**7*, 


it follows that its solution is obviously T**T,+(.T,-T,)e t 
dT 

Moreover, U ^r = (r.-T,)e‘ . If we take this expression and 

dl 

substitute it into equation ( 15 ), we then get 


(w).— 


(16) 


e substitute equation (13) and equation (10), equation (10) 


changes into 


djT dT £../ T,\r ( C+E\T._ FT. P 

df ~dl ST ~ R'\ l T, ) LV C IT, CE„ Pf 
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Its boundary conditions are: 


£= — ootime 7 1 = 0 
time T — T, 


Therefore a general solution for the equation is: 


T 


c'e 




T, R*T. 
T. CE t , 



In the equations R'T.V T.) • Moreover, is then 

determined from the characteristic equation A i_ A _s = o . Because, 

when £— oo , T = 0 ; therefore, one only has + + )/2 

as a characteristic limit. Again, from the boundary condition £=o 
, when T = t,, , it is possible to obtain a solution to the equation as 

follows: 


T 



T. R'T. 
f. CE M 



( e U_ e t) 1 

-J 


Substituting in equation (16), one then obtains; 


i+ d( c+ e.._T-_ E1L. t\ rA n 

* + S\— c - T. CE„ P) 


(5h 


f 1 71 











Prom equation ( 13 ), one has 



c+ e7. 


c-+i- 

ffi P 


( 18 ) 


If one takes equations (14) and (18) and substitutes them with equation 
(17), one then obtains: 



C+- 


’ •*_ 

P 


m P 




C + E., f m~ n 7 . r_i 

c™ 4.L i £/_ t, 
m P 2 R*f f 


C ^ + ? 
R'T, 
E.. 


1 


n- 2KT, 

If we make ~ (T.-T^E/ . and we make' use of the pressure response 

function definition R r = jL^ , then, the above equation can be 

arranged so that: 


n AB— 


A — 1 




C + q/E„ 


r-y4 + /IS+ —— . p. - 

/ C + q/E„ 


(A-l) 


Area II: 

Based on a method similar to the one above, we only take basic 
equation (2) and change it to be m = ^„exp(- -fiif-) • It is then 

p :m j’.ble ho )b t a i n an expression for the ores nice response function for 
Area 11(5) : 










/? n = —--- 

A-A + AB+ (A-l) 


Again, on the basis-of the normal assumption that there is a direc 
proportion between combustion surface gain and combustion surface area, 
the pressure response function for the entire combustion surface formed 
by Area I and Area II should be: 


R,=y- 


A— 1 

n 'A'B<-c + q/ E;, 

C 


T -A,+A,B,+ c + q/E t . 


(A-l) 


+ (l-y) 


it\\AnB\ 


A' -An + A, t Bu + (X- 1) 


Actually, the imaginary section is: 


?„n, A i Bi D + i A, Bi — c + 

(R,)r = y D’ + E * * 

. s n M ^iiBn(A,D + A|F) 

+ u v) D'+F 1 

— X r n, A t B ,EA, B i — C+q/E^ ^ 
U?#>).=*v D' + E' ^ 








D = A + XAAB-A), E = 


^_p + x^ab-a). 

F sQ + UAB-Ah 

A '-Wl[- ,+<i+i6 °' )i j i 


4 Results of Calculations and Discussion 

Waking use of the steady state model(l), the various results of 
intitial calculations on S04-5A experiments (Table l) form the Initial 
data to substitute into the formula for the real section of the pressure 
response function, yielding results such as those shown in Pig. 3- 



^ig. 3 The Real Portion of Two Low Pressure Response Functions for 
S04-5 A 


The formulas obtained and the results of calculations all clearly 
show that the real sections of pressure response functions for "plateau" 
and "mesa" propellants whose pressure exponents are aero or negative are 
all capable of being positive values larger than zero within quite a 
large range of frequencies. From the standpoint of theory, this explains 
very veil facts associate! 4 ith severe phenomena of unstable combustion 
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and the existence of these two types of propellents(4)(6^(3). Moreover, 
just as this author calculated in reference (1), due to the fact that the 
abnormal combustion of oxidizer covered with molten binder agent takes 
the place of partial shut down, in this way, under the pressure 
associated with the existence of widespread covering, it is, on the 


contrary, e 

oscillation 

asy to 

. It 

see the appearance 

is then possible to 

of the results of s 

give a reasonable 

elf-e xcite 

explanation 

* 1st 

<r. 

’ «c * ffi 

*1 ft 

i 

* a fit 

n 

; JTffl 2 

1 

B, ; 

kg/cm’ 

3. 48 X 10 s 

r 

kg/cm’ 

23 

1 35 

c 

Cal/g c K 

0.3 

G 


0.3 

0.2 

E a , 

C«1 mnl* 

2.2 X10* 

V 


0.6 

0.8 

Ei 

C»1 mol* 

j 6.0X10* 

- 


0.572 

-0 693 

fl i 

Cal. rr.ole 

2.08x10* 

T ,i 

•K 

1038.6 

1098.8 

r, 

" K 

! 2900 

*>/ 


0.75 

-0.783 

T, 

‘K 

293 

T,„ 

°K 

1086.6 

1108.0 


g mol* 117.5 :. n,, 


0.486 1 0.424 

HI g'mol* 1 28.4 1 



i 

f Cil/cm sic °K 1 0.3x10"' j| 

1 

f g cm' 1.665 

1 


Table 1 1. Parameter 2. Unit 3- Numerical Value 4. Sample 

1 simulation 

Vhat the dotted line in Pig. 3 shows is that, in reference (1), S04-5A 
propellant, with a large surface area covered and a pressure of / J = 23kg/ cm’ 
, in experiments with the "T" type engine, shows the appearance of 
self-excited oscillation. The frequencies involved correspond to the 
real sect ion values of the r ■> eponsa functions for q u 11 e high pressures. 

^ -»ca is— of this, the appearance of self-excited oscillation in these 
'■‘Xperiments cannot be considered surprising. 
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Most of the results discussed above were obtained through analysis 
in 4rea I. This was done through the use of our steady state model and 
with consideration given to the covering of the oxidizer surface with 
molten binier agent and to the existence of a congealed phase reaction as 
well as opposed gasification. Because of this, the theoretical formula 
derived for the pressure response function can not only overcome the 
deficiency associated with the fact that the pressure response function 
formulas we already have cannot be used with negative propellent pressure 
exponents, but also, in one respect, it reflects the vitality of the 
stable state combustion model(l) for solid state propellents, which 
presents the combustion speed characteristics for positive and negative 


pressure exponents. 


The formula which is derived when the covered surface area fraction 


is zero, then becomes the "two parameter formula" for 
propellents with binder agents which are difficult to melt: 


A + — — (1 + A) + AB 


lading .is e of the unified expression we already have for the "gas phase 
quasi-steady state, uniform propellent unidimensional model"(2), this 
article consequently derives a formula which is capable of having even 
broader applicability. 


This article, before deriving this formula, made no small number of 
empirically quite significant assumptions- for example, the assumption 
that the congealed phase layer which gives rise to the congealed phase 
^action and the gasification reaction is an infinitely thin plane is 
certainly not a good approximation^). The assumption that tne gas phase 
a 33 ua ?s a steady state has an error less than i 0% only when the 








j 













fre-.iuencies are smaller than 10,000 Hertz(2)(5). The compressiblity of 
propellents also has an influence on response functions which is within 
10#(9)* Heat radiation losses, under certain conditions, also have an 
obvious influence on response functions(10)(11). As far as congealed 
phase reactions are concerned, in our stable state model, they are 
basically pressure functions(1). However, we also conveniently saw from 
our mathematical treatment, that, when their influence on pressure 
diminishes , error is, naturally, even more unavoidable. Although, 
except for congealed phase reaction problems, the majority of assumptions 
are obvious, they are also universal in their application. However, this 
explains the fact that the formulas derived above are still awaiting 
refinement. 

At the same time, formulas derived on the assumption of homogeneity 
can, of course, only be used in situations where one has homogeneity. 
However, due to the fact that propellant impurities have already been (57 

used by several scholars (12) in methods such as the one that follows, 
that is, taking a multi-mode coupesite propellant and viewing it as a 
finely dispersed oxydizer granule/fuel "surface match up" of random 
arrangement. When one assumes that each "pairing" is mutually 
indpendent, then it also bscones possible to take the propellant surface 
and rearrange it into a family of hypothetically dispersed units of 
propellant- This is what is called a hypothetical propellant. And, the 
pressure responses of these hypothetical fuels are calculated using 
formulas deduced from homogeneous theory. Moreover, from this one 
obtains the pressure response for the whole composite propellant. Just 
as is pointed out in reference (3), although this method has the 
drawbacks of two different models in its application, it is still 
currently being followed. Because of this, this model certainly does not 
fail to explain pressure pairing problem values for composite solid 
propellents because of its adoption of a homogeneity hypothesis. 

5, Conclusions 

I:.:' ij i; - > n t 5 S Judy state coabustion uodel for positive and 

negative pressure exponent combustion speed characteristics for solid 
propellants, we deduced a new pressure response functlia formula. This 
^ormala explains, in a reasonable way, phenomena associated with unstable 
nimbustion which exists in all solid propellants with positive, negative 
and zero pressure exponents. 
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2) In the sHi 11 3 way which was the case with the forming of a model 
for steady state combustion, in the combustion of propellants, the 
phenomenon of molten binder agent covering oxidizer surface area should 
also be considered in the forming of a model for unstable combustion of 
solid propellants. 

3) Using the unusually regular combustion of the congealed phase 
reaction and opposed gasification covered by molten binder agent to 
replace localized shut down is a precise and necessary method. This is 
true not only for steady states, but for non-steady states as well. 

4) The model presented in this article is still awaiting further 
experimental testing and perfecting. 
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SUMMARY 


In this article, under simulated engine combustion chamber p.'es-sure 
conditions, we carry out an analysis of and calculations for the dynauos 
of the omni-axial rotation of flexible jet nozzles. The important 
contents include: analytical calculations of measured data; analytical 
calculations for the center of instantaneous oscillation of jet nozzles; 
and, calculations of the moments of force of oscillation. Besides this, 
we also calculated the paraaeters that follow: the angle of oscillation 
and angular velocity of jet nozzles; azimuth angles of oscillation; the 
length of operating tubes; eccentricity of thrust; and, axial and radial 
displacement of jet nozzles under differing pressures of contents, etc. 

This method of calculation, after having appropriately resolved the 
p”obiea of installing telemetry pickups, is still suitable for use with 
hot test beds and is also suitable for calculations and measurements 
under conditions of static oscillation of flexible jet nozzles. 

1 , Calculations of Adjustments in the Amounts of Uisplaoenent and Angles 
of Oscillation 


(1) Selection of Measuring Equipment and Coordinate System 
The apparatus for measuring simulated axial oscillation of flexible 
jet nozzles is as shown (?ig. 1). The flexible connector head connects 

" e-sp so ti vely with the end ring and the lower flat plate. They and the 

container shell body together fora a sealed high pressure vessel. To the 
hot to a flat plate is firmly attached a rod (simulating the jet nuzzle). 

Its axis line and the axis line of the flexible connector head are 

congruent. Its point of contact on the Z plate is the origin of the 


content pressure 


zero. 


The axis OZ is then congruent with the axis line of the rod. 

The top end of the rod and the axis line are installed at 
•’' gh t tngl e s f •> th * "1 a t plate, wh I. y\ ' s called the Z plate. The axial 
line >f t’-i t 1 i splac e meat sensor z, is congruent wita the axis line of 
th e od . "”i => direction end amount of thrust are coot"Oll-d by actuator 
tubes ' and 2 , which are r > spec tl vely positioned in plane 3 Ox z 

and OYZ • I n this way a fixed coordinate systei O .Y YZ 








is then established. Two flat plates are installed on the rod. 

They are respectively perpendicular to the X axis (called the X plate; 
and the Y axis (called the Y plate). Displace meat sensors 
and Yj.K, are respectively perpendicular to the X 

plat- and the ^ plate. Moreover, lines extending trot A. 

and Y, , and A', and V, connect respectively with two points on 

the 7, axis. Their points of intersection, with coordinates in these 
dimensions, are (0,0,-Z,),(0,0,-Z,) . An extension of tn j 

axis line of the displacement sensor intersects with the OZi 
coordinate plane at point (0, — Y,, — Z ,) . The distance from the 

plate and the Y plate to the coordinate point is F. in both 
cases. The other two displacement sensors H>- and H r.-jsp ecti vely 
measure the two actuator tubes fixed to rod points Q, arm 0: in 

terns of the anoint of axial displacement caused by elastic deformations 
of tne vessel. 

Besides this, the radius of measurement contact heads for 
displacement sensors X,*X , is F, in all cases. The 

radius of the measurement contact head for sensor Z, may be selected a 
will. Let this value be ft, . The rest of the ins ablation constants 
and structural constants are a; shown in the Figure. 

F, (or p, and F, (or ) are activating forces (or 

pressure differential) sensors. 

(?) Analysis of Displacement 3ensor Measurement Data 

Amo ac t Vi a eight displacement sensors, except for the measurement 
•1 f r > ! H, n’ ff 2 1 tich there is no need to correct, the measured 

.lit i )h ti i e d from the other six sensors certainly does not perfectly 
reflect the aaount of displacement in the oscillation of the rod. This 
is due to the reasons set out below. 

(1) The measuring contact head of the sensors is a hemisphere. 

A:ling t) changes in the angle of oscillation, contact points 
experi-Jic 1 displacement along the hemisphere, introducing an additional 
amount of displacement. 
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Fig. 1 Schematic Diagram of Test Equipment 1. Plate 2. Rod 3. Lower 
Flat Plate 4» Upper Compression Ring 5* Flexible Contact Head 6. High 
Pressure Container Shell 7. Actuator Tube 2 for Controlling the 
Direction and Amount of Thrust 8. Plate 9- Diagram 10. Plate 11. 
Actuator Tube 1 for Controlling the Direction and Amount of Thrust 

(2) The measurement contact heads for the various displacement 

sensors X,,X^X,>Y,*Y, and the contact points of the 

experimental me a.«nr ? .-jut plates are at a distance from contact points on 
the rod axis line X , and from the symmetrical surface of the rod so that 
the distance varies with changes in the angle of oscillation. 

(3) Because of the fact that the axis lines of X, Y, and Z rotate 
at the same time, their mutual influences cause plates X and Y to 
correspond to the slant which occurs in the sensors. This introduces an 
amount of non-displaceaent in th i sensors. 











(3) Adjustment Calculations for Angle of Oscillation and Amount of 
Displacement 

Let the angles of rotation for the axis lines X, Y, and Z be 
. In order to simplify analysis and calculations, we maK:e use of 
the geometrical principle of superposition to make a two step 
calculation. After that, we superimpose. That is, first we make 
calculations for the situation in which a= *0 and (J*0 , withfl«*0 

Then, we make the calculations for the situation in which 
a * o,/3*0 and o 

(1) Adjustment Calculations for the Situation in Which <*=*0,0*0 
and we assume 0=0 

There are four types of factors in omni-axial oscillation which 
cause adjustments to the displacements measured by the displacement 
sensors: (a)a>0,/J<0»(b)a<0, /?<Oj(c)a<O,0>O|(d )a>0,/9>0. 

To'’■aula deri/at Ion brings up one type of situation. 

Eccentric rotation of the axis line within the first apparent limit 
a>O,0<O. 

If we assume that the cross section of a rod is a square 2R on a 
side, then, Fig. 2 is a cubic diagram of angle a i\A angle # of the 
deflection. Fig- 3 is a horizontal projection diagram corresponding to 


Definition: The plane which holds the displacement sensors X t 

and Y, (or X , and V, ) is the horizontal plane of test 
measurements. The same rod cross section used in it is reprsented by uEFGH 
If we ussu'ie that axis line Y rotates through angle a and 

axis line X then rotates through angle P , the change in the horizontal (^1 
plane of test iieasur ements is: squar° DE,F,G,H t ~ elongated square 
zsE'F'G'H quadrilateral zjEFGH . This corresponds to a 
change in the contact point of the experimental measurement head of the X 
displacement sensor: A,—A,-*A lt 

At point a„ • a = 0,0 = 0 - X sensor measurements 

give us a displacement amount X, = o ,A,0 = R 

At point A,. a>0,0 = 0 , The nsccrb amount 

disolucement is X'. The actual amount o f d Isplao e -lent Is X ' = 0,0 


The mea'u” j d amount o J - 


One of the additional amounts of d i splace a-^ni 


A l O l -A,0=R/oosa-R 












Pro (i the diagrams it is possible to see four amounts of additional 
displacement caused by the displacement of the contact points on the 
contact heads: A C A, =0'A, -O'A^Rt/cosa —R* 

X' = X + (R+R,Xl/ooea-l) (1) 

At point. Ai. : a>O,0 <o . Let the numerical data measured 

by the sensor be X" . -From Fig. 2 or Fig. 3 it is easy to see that 
the value of the difference between X" and X' is A','A,, 

These are three of the amounts of additional displacement caused by angle /?. 
on sensor X . From Fig. 2 we obtain the set of relationships 
set out below; 

From RtAA x A\A\ we get: ' A',A‘,' — A,A',tg0 = Ytg/3 

From RtAA,A’,A t we. get: A'lA 1 = A' 1 A'itga*=Yigatg/3 (2) 

Froa RtAA,A,A', we get: A l A 2 =A,A t tgr — Ytgr ( 5 ) 

Compared with equation (2), equation (3) has: r = arctg[tga • tg0] (4) 

Because A' 2 A 2 = Yigr<0 > W e also give consideration to 
satisfying the geometrically equivalent relationship, and we then have: 

' X* — A'tA t = X' (5) 

After we now do some more analysis of the deflection angle 0. , 

the four additional amounts of displacement which are induced in the 
displacements of the contact points on the measuring contact head of 
displacement sensor X can be seen from Fig. 4. The contact point moves 
f -om point A,. to point A ,, : that is, it rotates through angle X.A'-0*A, = r,0*A l 

= R,cosa, Because of this, there is always a positive 
mount 0 " increase: 

J?,oosa(tyeosr ,l M) ' (6) 

When one combines this with equations (1)(2)( 5 )(4)(6) one gets: 

X-X'-Y V r-W + R.^- l )-R.c*a(-i--l.) ( 7 , 

Below, when we 'flake adjustment calculations for the amount of 
d Is pi ..cement measured by the Y direction displacement sensor. We 
similarly have: 



At point 



At point 1 
At point Bi 

■ X'tg r <0 


F=0 5.0 = /? 

y/_K+(/?+/?.>( i _j\ 

'«*/» ) ( 8 ) 

B',B,^Xtgr 
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it satisfies the equivalent relationship 
Y" + B'iBi—Y' 

Y”-Xigr = Y' ( 9 ) 


The additional amount of displacenent caused by the angle a is: 


/?,cos/?( l/cosa—1) (10) 

From the various equations (3) (9) (10) one then obtains: 

Y^Y* — X tg r —(/?+/?,)(1/cosl) — /?,cos/9( 1 /cost — 1 ) (11) 

(2) Calculations to Adjust the Amount of Displacement When the 
Angles Are Simultaneously Wot Zero 

As far as an ideally flexible contact head is conerned, one has 
only the times when the azimuth angle of oscillation is 45 degrees 
and 225 degrees , making use of forces, the combined force of which 

crosses over the axis of the jet tube. In the case of oscillations with 

other azimuth angles, making use of forces the combined force of which 

not at all cross the axis of the jet tube,adds an additional moment 

of force to the axis line of the jet tube. In this way,the jet tube 

rotates around the Z axis line, the angle of which is designated e 
During omni-axtal oscillation, it is only possible for the angle 9 to 
be positive or negative. If one gives consideration to a,/9 , the 

overall influencing factors can be grouped into six types of situations, 
that is, 


0° *= <45° 

45°<<t.<90' 
90° < ^,=£180° 
180° < V,^225° 
225°<^,sS270° 
270“<v>,sS360' 


0<O,(a>O,5<O) 

0>O,(a>O,/9<O) 

0>O.(a<O,0<O) 

0>O,(a<O,5>O) 

0<O,(a<O,5>O) 

0<O,(a>O,5>O) 









'S'ig. 4 . Movement of Contact Point on Measuring Contact Head 

Pig. 5 Diagram of Rotation of Horizontal Test Measuring Surface «<o 

Doing through a practical derivation clearly shows that the six 
situations all possess the same type of expression. In this case,it is 
only necessary to do calculations for one situation. 

O°^^,<45',a>O,/3<O,0<O,A'>O,K>O 

0 is the angle of rotation around the OZ coordinate axis. It is 
) 11 Vi it th ■> h 0 •sental test aeasurement surface rotates through 
ingle 9 in the same way. Analysis of the front surface down to the 
mutual influences of the angles of deflection a and 0 cause the 
horizontal test measurement surface to give rise to deform it ion. After 
maxing adjustments, a and 0 are used independently to cause the 
horizontal test measurement surface to give rise to deformations in o E'F'CH* 
as in the diagram, and rotate it through the angle 0 • 

th >t is, becoming zdE"'F'"G'"H'" , as shown in Pig. 5- 

?roa s'ig. 5 one can see that the angle 6 causes the X sensor and 
v 1 sn s >r t > r •» so a ctivel'r aenui re an additional amount of displacement 

V A. X 

beooiing A','A, and 
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Solving for A'iAu , the D point coordinates are 

{X o—X +O t DccsO 
W d **Y — O,Dsm0 

A.O.-O.D,, 3 LA\0 = X' = X+(R+R.)(— -iWa'+^I.O, 

\ cos a / 


:.0,D=A,0, = (R+R.)( — I- -l) 

\ cos a / 
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, also because 
therefore, the D point coordinate should be: 


DC~0,D sin0<O , 


' (/? + /?,)(—--l\cos0 

I V cos a / 


( 12 ) 


Make 


and H'E' cross at point J. Then, the equation 


For the straight line DJ is: 


Y —$i <-X-X.)+Y, 


hi. so one iiae: 


X,= X' = X + (R+R ,)(—i-1) 

\ cos a / 


(13) 


(U) 


Take equation (14) uni substitute in equation (13): 

y,--±\ (R + R,)( — —-l)(l-cos0V + Y + (R + R,)( --ljsini 

\ cos a / j \ cos a / 


( 15 ) 


^rO'Ti Rt.\A':Aj we get: 


A ; A, = YjigO 




( 16 ) 
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Prom Pig. 5 we get the geometrical relationship 


a,o+a‘;a 3 =a;o=x' 

A,0=X' t A';A, = y, tg0<O :.a,o-a:a,-a;a : =a,o 


Substituting the various algebraic quantities and simplifying, we get: 


X = X’-Y(tgr+tg0}-(R + R.)( 1 _!) 1 _*,«*«( 

\cosa /cos# leas r / 


(17) 


In solving " >r B'jB, , we use a process identical to that use! in 
solving for A"A, , and we obtain an 3 point coordinate which is: 


[*»-*-<* + Mjshr ') sln9 


(18) 


Let E"'F'" and £"F' cross at point Q. The straight line equation 
of 3Q is 

K*t g fl(x-jf 5 )+y. ( 19 ) 

The Q point coordinate is: (20) 

Y'-Y-Y^R+R,^- ij-l) 

T it ? illations (20) (id) and substitute in equation ( 19 ): 


o a 


RtAB.B'iQ 


L 3 possible to obtain B ",B,= . Yotg# 


( 22 ) 


The geometrical relationship has the form: B,0-B\B,= B-0 
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■.B‘:B, = XoigO<0 


; at the same time, after giving consideration to 
the superposi tion quantities for B'iB, , one then has: 




(23) 


Equations (17) and (23) are the equations for calculating the 
amount of displacement adjustment which is measured by the displacement 
sensors in the X and Y directions during omni-axial oscillation. The two 
are both functions, and, set up together, they can be solved for an 
equation through which the amount of displacement (X*Y) can be 
calculated from the measured amount (JY * % Y r ) : 


Xn = 


ccsfl,\ cosa. ) \cos0i ! 


(tgr, + tg^«)J l)[cosa 8 -cos^ ( (tgr ( +tgfl,)] 


(24) 


+ ig’9i 


Y',:-x::(tgr i - i go,)-(N+N t ) * r ( * ±~-i) 

y CO60, . \ cosjS, ' \ cosa. ' 


1- ig'r, 


(25) 


(tgr. - tg0,)] -/?, ( -— 1 )[oos0,-oosa.(tgr,- tg0.)] 


+ tg '6i 

!n the equations ;=l,2,i —0,1- 




are calculated from the 


Ah/?o- constants; 

formulas given below. 

(3) Calculation of Angle of Oscillation 

From the definitions of the various angles of rotation one can 

have: 

tH a x ,. = x;:-x;;-(K ; ;-K;:)(tg a,tg^. + t g 0,) 

Z t -Z , (1 — tg’a.tg^i + tg'fljz, (26) 


,„/?=_ A'?. rJ'_L = YV'-Y': + (x\:-x',’,)(ig a,ig 0 ,-ig 9 ,) 
Z t Z, U — ig'a.tg'/S. + tg’fl.OZo 


(27) 


Calculation of the Angle 


o , 


From equation (17) it is 


possible to obtain: 













TWTW7* WW Iff Iff IF ! IffIff WHmnMW1I I IHiW»W fflff iff w u i ipwnwmwwww 


-^ 2 i *X;:-K li (tgr l + tg0 i )-(/? + /?.)f—--1)— L-A'.cosai-^-l) 

V CXKfli / OOs9i V COSTj ' 


Through deduction it is possible to obtain: 


A r J .=^;:-(K l< +K.)(tgr i +i g 0 1 )-c/?+/? i )(^—iV-- 4 - ' 

\ cos a, / cosfl, 

— /?,cosa/ — -l) (23) 

V oosr, ) 


? r o i tbs 1 e n. a i. t i o 'l: 


, _x>,-x„_x',:-x;\ i._. 
18 —y~- 2 ?:— 2 ,gr ‘ 


(29) 


liquations (26)(27)(29) are a set of implicit transcendental 
equations. Using the iterative substitution method in a computer 
application, one solves for values of a<.£,-.0j . After this, one 

then substitutes equations (24) and (25)- It is then possible to solve 
for X,i and K,< for the various amounts of displacement. 

(4) Adjustment of the Numerical Data Measured By the Displacement 
Sensor 2, 

Definition: The angle of oscillation is the angle included between 

the geometrical center line of the jet tube and the OZ coordinate axis. 

It is expressed by the use of i 

The numerical data recorded by the 2, sensor, alone, due to the 
-ngl -» of slant on the Z plate, introduces an additional amount of 
Uspinesaent. Moreover, this angle of slant is always equal to 
oscillation ■ ngls $ . In order to solve for this oscillation angle, 

ois oust ?L '.}t solve for the direction number of the geometrical center 

1 i_ Q 3 


A> = A ,i — X ,i i Bt — Y u — Y,i t C = Zj — 2, = 2 0 


(30) 
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Cb5 











Oscillation angle t 


coe o ,= 


*J~A\+B\+C> 


(31 ) 


Adjustment calculation for ZJI ; Z u — Z\\ — R ,(— -(32) 

\ coetf. / 

In a heat test bed situation, it is possible to make use of two 
d i 3 pi ace nen t sensors Z%a - uni Z,» la order to replace sensor Z» . •< 
uni Z,, are installed symmetrically on the two sides of the jet tube in 
order to avoid ignition flame heat corrosion. 

As "Pig. 6 shows 



' IS 
a ni 


v a, 


Ohart of Sensor Distribution (b) A Comparison of Sensors 


to Sensor z. 


* \A 


Z, 


and Z, 
(33) 


It is easy to demonstrate that, between 
there is the relationship shown below: 

z;'=(z;;+z;;)/2 

(5) Adjustment of the Coordinates of Actuator Tube Piling Points 
Cj'it'ol the Direction and Amount of Thrust 

Due to the results of pressure contents and the forces used with 
actuator tubes, the displacement of the rod filing po : nt >.lotj the at is 
reaches above 5m i. This not only produces an influence on tie shock 
process of the actuator t foes. It is also related to the ii-eotiou o f t 
forces used and the force momentum. Therefore, it is 

necessary to make adjustments to the measured amounts. Host of the 







amount of the oh a nge is along the direction of the axis. Therefore, the 
adjustment is approximately: 


(Zqii — H q + H\ i 

^ Zqh = H o + Hit 


(34) 


In the equations, Hq is the coordinate of the rod fixing point. 
a re the values of the amounts measured. 

(S) Calculations of the Motive Forces F, ( and Fu Used to 

Control the Amount and Direction of Thrust F, 

As motive force one can choose to use the tensile pressure sensor 
contact measurements obtained as F, t and F,i . Or, one can measure the 
pressure differential P on the two sides of the piston. The force used 


is then: 


Fn = SPi 


In the equation, S is the area of the piston, 1,2,- n. 

Up to the present time, we have solved for the following parameters. 

dr, t*-Xti■*-Xit\ ^ Za,i\Zotiyt 1 tt^Fn%oi \ ftt 1 \ f t \ it yAt\lit 

and C , as well as constants to be used in later calculations. It is 
not necessary to explain these for our purposes here. 

2, Calculation of Center of Oscillation for Flexible Jet Tube 


Tl-C CTiOdPT OP A CCUTUI )P 030 [ 60AT [04 

As far as the concept of a center of oscillation for flexible jet 
tubes is concerned, both inside and outside China, there have been 
treatments in a amber of technical reports. But, these reports have not 
found a clear and precise definition. Moreover, there have been several 
methods put forward which are obscure and unclear. Before making 
c ilealat ions of the center of oscillation, we plan to have a clear and 


precise concept to act as a theorm. 
calculations. 


a i -ta . i- 


or analysis and 


■m t j i'\ ^r 'j i 


uder stood in the 


following terms. In the case p* rig' i md'.es, there are, of course plane 
oscillations taking place as wo 11 a’. \ it', al m'illntions (omni-axial 
p sc illation) • Both of trvuse exist simply as a se r • as of instantaneous 
axial rotations, and no instantmmm center of oscillation exists. 
Obviously, a jet tube is a rigid body, and its center of oscillation can 
chanqe. 
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Flexible jet tubes also have a point which functions as their 
theoretical center of oscillation. If the design and nanufacture of the 
flexible connector heads are both absolutely ideal, this leads to the 
flexible jet tube rotating around the same fixed point. What are called 
its instantaneous axes of rotation all pass through its theoretical 
center of oscillation. 

Absolute idealization of flexible jet tubes does not exist. In 
reality, due to content pressure (or combustion chamber pressure) 
changes, flexible jet tubes are made to produce relatively large amounts 
of axial instability. Besides this, flexible connector heads do not have 
uniform distribution of material around their circumference, and, during 
oscillation, their acceptance of forces will not be entirely symmetrical. 

It is also impossible that the deformations in the flexible connector 
heads, which a- e caused by this, should be symmetrical and uniform. This 
will also cause the instantaneous axis of rotation of flexible jet tubes 
not to pass through their theoretical center of oscillation. When 
flexible jet hub as undergo omni-ax Lai oscillation under various types of 
content pressures, the series of instantaneous axes of rotation will 
form a cylinder-shaped included figure. This figure will be symmetrical 
around the center of oscillation which will be included in it. This is 
called the "oscillation center" for short. This is also nothing else 
except the actual oscillation center of the jet tube, the size of the 
geometrical dimensions of which reflect the qual ity characteristics and 
rigidity of the flexible connector heads. 

If one is considering directly solwing for the geometrical range of 
the included center of oscillation, it is necessary to have dense 
instantaneous axes of rotation. It is also necessary to have large 
amounts of experimental data. It is difficult to get to this point, and 
it is also not economical. On the basis of this special characteristic of 
jet t ib es -o tat lag around these fixed points, it is possible to reach the 
conclusion which follows. The instantaneous oscillation center of the 
gecaetri cal center line of jet tubes, is always a reduplication of this 
fixe! point. This p iX' j d point : ■? what is called the oscillation center 


of this, we then use 


instantaneous 


enter 


geometrical center line of 


jet tube to 


replace the point center of oscillation of the jet tube. By means of 







ijvyw'i* wt-'v: >*■■■» r* - - *v»vw. ^ i,’*',T\'»v'»i.v L >vwflr^JVWWP^ 


this method, we use a limited amount of experimental data to solve for 
the geometrical size range of the jet tube oscillation center. Moreover, 
this supplies a reference point for solving for the Lnstantaneous moment 
of force of oscillation. 

During oscillation, if the geometrical center lines cross each 
other, the results of the calculations discussed above are co mletely 
accurate. However, it is also possible for them to be mutually off and 
not cross. In such a case, one then taxes their "average oscillation 
plane" and solves for the oscillation center. Obviously, the results of 
this are approximate. The degree of this approximation varies with and 
increases as the degree of graduation of the measured angle o f 
oscillation is reduced. Experimentation cleanly sow-- that, vnsn the 
measured angle of oscillation is change! f-on '5 iegre-rs to 1 .5 degrees, 
the distance between two adjacent geometric! • ?uter lines drops from 
0 .5-nm to 0.2mm,roughly spea'xing. 

CALCULATION OP THE OSCILLATION CENTER OP PuEIiti.E -IE 1 TUBES 

In order to simplify the discussion, let Ye- gsonetrleal ; enter 
line of the jet tube be the Z' axis. This Is riot the saae as the 
axis of the coordinate system. 


The 


axis equation is: 


X-X„_Y-Y 11 _Z-Z 


A, 


B, 


(35) 


’he equation for the Z plate plane is A,X + B,Y + C(Z — Z,,) = 0 


The set of simultaneous equations (3o) (3b) gives us 
coordinates for the 0, point. The 0. p'int -Ye it ? 


, —CZ, l — A,Qt — B i P, 

‘Oj — 


A , A i -+- p 11 , + c 


y»i — L * • ( i 
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5 o ! 


v * s - v - * - y Ua V. vvV> 




Lv.;-', 
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In the equations 


K,^Ai/CfQ, = X,i — K, i Z,}Li = B i /CiPt^Y,, — l i Z l 


The 0< point is a point on the Z' axis. In the initial 
configuration, it is a duplication of the origin point of the coordinate 


system. 


If we postulate two mutually adjacent instants, the plane at 0 i which 


contains the instantaneous center of oscillation is x„ 


We then 


nave: 


(jc —x,<)’+ (> — >,. )’ + (^—*./)’ = ( x— jc„_,)’ + (v 
After rearrangement, we get: *u« A,iX+B,iy + C,iz=D, ( 


In the equations 


(Aii — Xti Xu— ifBn — ya )C,. — 2 tl 

<D ll = (x,J + y,| +z,\ 


In order to solve for the instantaneous axis of rotation of jet 
tubes, one must still find a second point on the plane on which we 
find the center of oscillation for the time interval during the same 
period. Let this be the intersection line of the plane *»/,*»< and the 
plane , and one then has the instantaneous axis of rotation. 

Leoause of the fact that a measurement error in the system is 
i.inavoilabl ^ it is only when »>.•!-» h. that the degree of precision 

associated with the axis of rotation which we find is optimum. The O, 
point is a Pixel point on the jet tube axis line. Points which satisfy 
the joiliti.oo will vary with changes in the oscillation of 

the jet tube. The method for the fixing of coefficients is none other 
than an unceasing pursuit of these points. 

If we take the plane t and set x« f n t , 


A,iX + B,,y + C„z — D m 


If we solve the simultaneous equations (35) all [}')) 
'j e: i 1 1 ' ^ ; j ; t i o i point of *«• and the Z\-, •<< 

When we solve, we get the result: 


E„(x„, y„,2„ 












’» _**,*« .*» ,*• 




•*> 















The coordinates for point £»• : 


in the equations 


' (D,i — A,,Qi-, — B _ is n • q 

2 “ AuKi-, + B,,Li-, + C,i K " *' 0,1 

x,, = K ,,D, + Q,, 

y», = A u D. + Qu 


Kt < L,_, + Cu) 

K K h K Lt-,K i, 

Oa “ Kt-iQf •r O.-i 0«. = L,-Q u + P, 


'A x'/i e the i jr> -5 -> jt ion point of plane and axis Z’, be 

Whin nl/3, ve ^it thi "isult: 


Coordinates for point ft* 


f : _ d..- a,,q,-b ,,p, 

" A,,K, + B,,L, + C,, 
x,, = L.D., + P., 
h,.*L...D.. + P., 


in the equations: 


A,,K, + B i, L, + C : i 

Li * = A (Liii 
ii, e Lfiii i 


n . n • 


M _ * B„P.) 

' " A,,K H..L,*C„ 

P„-Q, + K.P„ 
P„~P, + L,P,. 


K..L..P..U. 


• ,• A •. 

•. * • *" t y , 




.vv.y.hhi-. A'vA'vt/.s % - 
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4fte- one trices e quations* (40) and (41) and substitutes them in an 
expansion arrangement, one has 


in the ui nat ions 


U tl D m l + V„D mt + W ti = 0 

U;-ltl +1,\ +1,1 -K»] -K,1 -K 
V h ~ 2t K ,,(*«,-i 0 n ) + K u(yn -i Qu) + K„(z tl -, — Q, ( ) 
Lu{x»t P i<) Ln(yti / > i,) Li,(z ti — P ,,)] 
’"‘(xn-P.ty + iytt-y„)' + (:„-P„y-(x„-,-Q u y 


)l/i.ng equation (42) we have: 


D„, .- -Y.,±*/r.!-4U„W u 
2 U tt 


(43) 


is then the fixed coefficient we were solving for. There 

‘o ^ts of solutions. However, cor responding to our actual problem, 

to- --elution j et in which z„ and *». are simultaneously smaller 

sere is aooii ruble. Therefore, we take the two solution sets of 
1 * D.n.< 

h w - j > ’. / i 1 ~ ) a t 1 1 jibstitute the a respectively into equations 
ani 4*). ?roa tie principles discussed above, we decide and make a 
:t' on. Because of this, after we make a decision, the solutions of 
i •; 1 r :h ./ • ;h ) ? J , then nuke it possible to solve for the 

n,., on which is located the simultaneous oscillation center for 


Moreover, 


*». JLw,., 


( X - y + ( y - y %t )> a _, it )* = ( , _ )- + ( y _ Vt y + ( , #1 )» 


A„x + B„y*-C,,y^ D„ 











in the equations 


(A tl = x, { — x tl , B„ = y tl — y tl , C„ =* z tl — z,, 
iD„ = (x,l + y.I +2,) — jc« f — y,;— r,|)/2 

The'line of intersection between planes -*u and is the 

instantaneous axis of rotation we are solving for. Below, we also solve 
for the plane of oscillation of the geometrical center line Z' axis of 
the jet tube, and we let this plane be 

During omni-axial oscillation of flexible jet tubes, the most 
general situation for the track of the Z' axis is a curved surface in 
space. If we take this curved surface and divide it into » partial 
curved surfaces, and, respectively use a plane to approximately replace 
their various oscillation planes an. , it is obvious that, when the 
dissected areas n are increased in number, the level of the 
approximation is raised. We select the number n so as to satisfy the 
requirement for precision which happens to exist. In order to solve for 
the oscillation plane , one first does the calculations set out 

below: 

(i) The direction number of the common perpendicular line 
between the jet tube axis line ZI_, and Z\ for two adjacent 

ins tants. 

Make the direction numbers of the axis lines Z\ and Z.', 
respectively N,, and N,, , then: 

N x {A,,B l ,C)**{A, l ,B,i,C u ) ( 45 ) 

in the equations 

Au = C(B t .,-B l ) t B u = C(A l -A t .,hC, i =A i .,B t -A i B 

(ii) Make plane #<->* contain axis ZI-. . And, make it 

par allel to N ,, Then, the equation for is 

x ZM— if y~~Y i1 1 2 ■“Z, 

ftf-ii A t -,, C =0 

Ant B nt On 




• • % • * - *. * . • . - .• , • •* 


« , w . .■ 







Or it is: 


( 46 ) 


A tl x + B t( y + C„ 2 -D, 




in the equations 


At, B B ( .,Cu~^tiC, Bu~Ai t C A,-iC u 
C % , = A ( - f B li -A tl B t . l , D tt -A tl Xu-> + B>' Y «-' + C '‘ Z ' 


(iii) Make the plane contain the axis Z' , and make it 

parallel to Ni. Then, the equation for R t is: 
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1 x X ,i,y—Y llt 2 — Z, 
R it ■ A |, Bn C 

I An , Bt t t C u 


— 0 


Or 


A t ,x + B t , y + C u z - D ti 


(47) 


in the equations 


An “5,C,, B,,C t 
C, t —A t B n AtiBit 


Bn**AuC AiC u 

D tl = A t ,X ,i + B,,Y H + C..Z, 


The intersection line between the planes Ri-i and R, are the 
mutually perpendicular lines Z'~, and the Z\ axis because of the 

fact that the direction number ", of the line of intersection is 
equal to the sector area (N„xN tl ) of the direction number of Z;_ 

and the Z\ axis. 

The various simultaneous equations (35),(46),(47) necessarily make it 
possible to solve for the intersection point of the mutually 
perpendicular lines of the Z',~, and Z\ axes . The reason for 
this is that the intersection line of and R t -is located on 

on the one hand, and,- therefore, crosses the Z',~ x axis, but, on the 

other hand is also located on the plane R, , and, therefore, must 
necessarily also cross Z\ , as is shown in Fig. 7- 

(iv) In solving for the intersection point Ei-i of the mutually 
perpendicular lines and the Z’,-, axis, take the equation for the axis 
line Z;_, and change It to be of the form of a parametric equation. 

From equation (35) we get: 

iX^At-J + X,,-, (48) 

j y^Bi-.t + Y,,., 

1 z-Ct + Z, 
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Pig. 7 The Mutually Perpendicular Lines of the Geometrical Center Lines 
for Two Adjacent Instants Curing Jet Tube Oscillation 

If we take equation (48) and substitute in equation (47), and, we ! 

also make t = K„ , we then have: 

K„ =*(Dn~At t X i,_i — BuYC, t Z + B ti B + C,,C) (49) 

If we take and make iterative substitutions into equation 

(48) , we then obtain coordinates for the intersection point of 1 1 

the mutually perpendicular lines: | 


F. i— I ^ t— I X 4 | "t X 1 — 1 

ye<-i ,B Bi- l Kt, + Y i,—, 
z t ,.,-CK <, + Z, 

(v) In solving for the intersection point E, 
perpendicular lines with the Z\ axis: 

Change the straight line equation of the axis 

form: 


(50) 


of the nutually 
Z, to a parametric 


I x-A,t + X„ 
j y = B,t + Y,, 

1 z-Ct + Z, 

151 
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If we take equation (51 ) and substitute in equation (46), and we 
then make t = K tl f have: 


is — u>« AhX i CgiZ 

*' AM< + £„£, + C„C 


(52) 


If we take and make iterative substitutions into equation 

(51), we get coordinates for the intersection point E, of the mutually 
perpendicular lines and the axis 25 ; 


E 


/I 


*ti ^ A t K g, + X lt 
B,K tt + Y „ 

= CKt, + 2, 


(53) 


It is >asy to see that the points and E, are certainly 

not the same geometrical points on the axis line of the jet tube. 

(vi) Coordinates for the center point £„ of the mutually 
perpendicular lines E,-,E, 

, +**<)/2 

Coordinates for the point £ #< 1 + y* ( )/2 (54) 

Z t*i “ ,)/2 


From equations (45) end (54), it is possible to obtain the plane of 
oscillation *».• fro a instant f,_, to instant t, of the geometrical 
center lines if the jit tub*. Taut i 3 , 


or /i, l x + B,,.y + C I ,*«D, l 

in the ->quations 


( 55) 


(70 


are t n e 


same as before. 











The system of simultaneous equations (33)/44),( 55) 


i-A\i x -t~ Bny~bCn2^Dii 
\AitX + Bny + CuZ* = Dii 

^ A, i x + B,,y + C, i z**D, i 


makes it possible to solve for the instantaneous center of oscillation of 
the geometrical center of a jet tube . This is 

also none other than a point on the center of oscillation of the jet 
tube. It is possible for this set of equations to be solved through 

the use of a computer. Therefore, we recognize that x ti ,y,,,z t , are 

already known. Later, we will make direct use of them. Fig. 8 is an 
illustrative diagram of the instantaneous axis of rotation of the jet 
tube and the instantaneous oscillation center of its geometrical center 
line. 
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3, Calculations of Moments of Force 

The systems of moments of force which we are solving for here 
point, under differing designated wave forms, to a selected series of ^ 
instants (or, one might say, angles of oscillation) to figure out the 
overall instantaneous moment of force. The focus of this article is on 
solving for the mo lent of force-angle of oscillation function 
relationship (A/ = /uJ)) . Prom this, it is possible to obtain the 

maximum and miniumum (algebraic values) overall moment of force M■•> 
and , as well as its ooment of force-angle of oscillation function 

relationship. Besides this, it is also possible, through the designation 
of different wave forms (normal sine waves, sawtooth waves, and square 
waves) to obtain the curve of the functional relationship M=/(0) . As 

far as the elaborate separating out of elastic moments of force, 
asymmetrical moments of force, and moments of friction, and so on, is 
»rice rued, this article makes no additional discussion. 

After one solves for the instantaneous center of oscillation, the 
major contradiction in a solution for the instantaneous moment of force 
becomes the problems of the point of action and the direction of a system 
of forces in a specified space. We chose to make use of the two methods 
of Euler transformations and four dimensional numerical transforms, and 
the relative difference in values of force arms which we obtain in the 
solutions is not greater than 0.8$. Even so, we only introduce one 
method here. 

(l) To i Print of Action ijsed to Control the Direction and Amount of 

Thrust 

From Fig. 1 one can see that the actuator tube movement points are 
~ispaoVvely G, and G, ,the points of action for motive forces 
and / r i . Moreover, the points G, and G, are connected to 
the rigidity provided by the rod. Therefore, the problM V>ic>us one 
solving for the spatial coordinates of points G, and G, 

(l) Selection of a dynamic coordinate system 
(a) Translation coordinate system: 













GX'Y'Z' 


-the coordinate axis which corresponds to the fixed 

coordinate system QXYZ is always maintained parallel. The origin point G 
is the base on the rod center line for point G, (or G, ). 



Fig. 9 The Transformation of a Dynamic Coordinate System From GX'Y z-gx'Y’Z> 


(b) Potary Coordinate System: 

GX'Y'Z' - the angles of rotation corresponding to the 

coordinate system GX'Y'Z' are and their rotational order is as 

Presented below: (See Fig. 9) 

First, we rotate through angle <p around axis GX' , obtaining 

coordinate systme GX,Y, Z, ■ 

We also rotate through angle 6 around axis GY, , 

obtaining the coordinate system GX,Y,Z, 

"Finally, we rotate through angle ip around axis GZ t » 

obtaining coordinate system GX'Y'Z' . The angles <p,<5,V> are 

called Euler angles. 

(2) Euler Ingles and Their Transformation 
(a) Solving for Euler 4ngles 

In the calculations of Part One, we already obtained the result 


that the angles of rotation around the various axes of the fixed 


coordinate system OXYZ by the axis line of the jet tube are 

’Tow we must solve for the Euler ingles which correspond to each instant. 

^roa ^ig. 9 it is possible to obtain the transform angular velocity 
projection relationship Tor the GX'Y'Z' coordinate system. 













cosd'tinp 


cos & co *W 


^ro«n Table 1 we get the following series of relating equations: 


dfi d<p , - 

a -—£‘ ir +s,nl -fr 


da di 5 . . (ftf 

6 j y '= -r- = coegs • - coeosini/i - -jf 

or df dt 


' d0 di 5 , dti 

co,• = -j 7 = sin gj • - 7 - + cos 6 • cos ® • —~- 

„ di at at 


! d 0 *® d<p -f sin 3 • d tp * 

da^casq) • dd—cosd • smq> ' df * ( 56 ) 

c/0* sings • c/d + cosd • cosg> -tty 

Coordinate transformations are only related to the results of 
rotational movements and do not consider the process of their rotational 
aovements. If one desired to solve for the corresponding angles in a 
certain situation, it is possible to recognize that one always begins the 


rotational movement from zero degrees arid simultaneously malt as us *. o 


method of freezing the coefficients, recognizing the use of approximately 
uniform rotational speed. In this way, equation (56) can be seen as 
becoming a set of differential equations with constant coefficients and a 
zero integration constant, for purposes of integration. Because of this, 


one has: 


, + sin d, • 

a ( = C06g?, • d, — cosd, sin<p, • 

0 , = sing:,d, 4- cosd, • cosip, • V. 


(57) 


In order to facilitate the carrying out of iterative substitution 
calculations, we change these to the form shown below: 


'■c 


















V( = (fl<cos?) < —a 4 • sin<p,)/cos3, 

<P. = £.-*< sin <5, 

<5. = (cr. -4- cos<5. ■ $in$\ • i»,)/cos <p, 


( 58 ) 


<p\ =a ( /coR95 lt 

The initial value for iterative substitution can be taken to be: 

/\ 

• finally, we solve for the Euler angles; & i9 

(h) Euler Transformation 

The purpose of the transformation is to carry out a solution for 
moments in a fixed coordinate system. On the basis of the transformation 
order which is given in Pig. 9, we obtain the results shown below: 

Transf ormation of GX r Y" Z 0 -*G X ; Y,Z, 



— sin i >, 0 

corv* 0 | 
/ 


(59) 


Transformation of 


GX,Y t Z t -~GX,Y>Z t 


cos 5 0 sin <5 

A, = \ 0 l 0 j 

— sin 5 0 cos <3 

Trans for nation of GX,Y,Z,^GX'Y'Z' 


( 60 ) 


Ai m 0 cos cp — sin v 
L 0 sin V' cos <p ■ 


61 ) 


be 


ti? tr visfo-’metion natrix from 
A ,, Then, we have 


GX'Y'Z’ 


to 


GX'Y'Z' 






I ccsiptoosd, , — sin ip,cxxd, , 

A, = , sin ^,oos<p, + cosv<*n d, - simp, , cost, * cos<p,— sin ^sin^sincp,, 
sin^, • sin^, —cosjt,sin(5 ( • cce<p ( , oos^,sin?>, + sin^, • sin3,cos<p (> 


sin <5, 

— cos 5, sin <p ( 
coe3(Cos?> ( J 


The transformation of point G,, in coordinate system GX'Y'Z' 


r x cu *» r &G i 
I yen \=A,‘ 0 | 
'* z'cu J L 0 J 


The transformation of point G„ in coordinate system GX'Y'Z' 


x cu ~i r M 

y'cu ~A, R c ! 

a„ J ' . J 


(j) Tr ins fo-’ nat Lon of acting rod points G, and G, In a fixed 
coordinate system 

As far as the translation operations of rigid holies are ? in 'em- 
the a in ant a* displacement of each point a i ?t b -* e vni vnl-n* . N ■> ■< , o 
i r 1 e ” t > sol^e for the amount of displacement of the origin pint c 

.)" the translation coordinate system: 


















Solution of the simultaneous equations (35) and (65) allows us to 
obtain coordinates for the point Q : 

^ZG.-i-V,-y/V*-4U t lV e ,t2U l ' —'^ 

s *<?,=■/f,*c,+Q, 

( _ T _ . D 


in the equations 


yci-LiZc.+P, 

U l = K',+L', + \>V i = 2lK i (Q.-x. i ) + L,(P t -y. i )-z. i l 
W i-(Q«-*.,)* + (/*i-y.,)’ + H l 


By putting equation .63) and equation (64) together respectively 
,-'.th equation (66), it is possible to obtain the t^ansfor nation of point 
(j, A*4 G u in rne fixed coordinate systea. 




(67) 


*c., “X CS) + x r , 

- V C!1 = , V Ct| +>■(., 
7 H = i Cl, + ?c, 


( 63 ) 


■ tion of the forces acting to control the amount and 


' . - 1 '-'hn 'sf forces acting 

no 1 ;p; i 'ate the axis line of the activator tube. 

, - ■ 11 nates for the dynamic fulcrum and the fixed 


h" their directional cosine. 

(69) 

4(i”Xcn RoiB h ~ yc ii= z ctl — H<)>,]> 
-^•i = *ci( • B ti = ycu Ro \G ti =*7c,, — hi Qu 1 






Arm length of activator tube: 


J.i-n/^.!+5.:+C.! 

+ B,! + C,! 


Directional cosine of forces acting: 


170) 


/• C 08 a,, = A ti /L 
F,,j cos^u = B t ,/L mi 
oosy. ( ™ C t< /L„ 


r cor a t , = A t ,/L y , 
als0 F„j CO 80 „ = 

— C tt /Ly , 


■; 711 


(2) Projections of the forces acting on the various coordinate 


axes: 


#F,„-F„C 08 a„ 
|F,»i*F,.-00801, 

F,„=F, ( co8y,, 


also 



F„cosa,. 

F„cos/?„ 

F„cosy„ 
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(3) Matrix Calculations 

(l) Solve for the matrix of the instantaneous oscillation :enter o 
the geometrical center line of the jet tube 

(a) Radius vector of the instantaneous oscillation center of the 
forces acting 

?u~ C,(? u and ~Ci&u or 


= (*6,<-*c, )/1000 

= (>'cu-0'r,)/1000 

also 

^ r i*i Jfri)/l000 

•j r »j', = (.y<?«, — vc,)/iooo 

= (zcu — zct)/ 1000 


r »«r = (- ? c»( *c.)/l000 


(b) Overall moment of force—Main moment and its relative ’o-iinb 
' * > j a jt ions of the main moment on the various coordinate axes : 







\ *. 


^ r r ,F,,,-',, F r - l r . , ' ,- r r F 

t 

M h ” ^ ( r >,ih - r „f ,, )•( r „/•' «. r 


' ..F , -r .,/• 




• v W. ■ . it , ’ ♦ W. 


\f, - vf. . #. 


i « 




A..-H ,C,.- r ,A„ - A,,C. “A, ,F,. - B.,A, 


ex* r_, 


_A/ ..A.. + M,,H_, + M.,C , 

^m.: -f a/7: i- a/,: % . 4 .: +5.: + cv: 

















jii* j-m ' l" 1 ) and (l~i / , we can then solve for the no neat of 


f *,n<* i nstuntaneous -ui j if "it it i 1 on 


M ., - ! M ,,! 


• i i * 


\ft ip-* •. f . u j i»eet of for oe : 




■ < . i i ' i 


; 'i e in ’ t• i - . i • A < t ) r h > t 11 i on 




i• - ■ ’ iti in “-nter which we already solved 

•i I , . t i ; t n e u possible to solv* 


■ ■ ■ •. i . i i • i i ; i < ; i 1 ' •* ) t a t L i n: 


_ 




, • ' j y o * n rough the horizontal p r o j - e t L on di^ran 
t < . 1 t r i*/t *•>n wni-h appears in the X— Y diagram, 

. ■ v ’ ~ \ i ■ : *nter if oscillation of the 


■i a’ , v * > ii * a i n its projection equation: 


■ — K „,x + b, 


K B„,/A„,.b, - y, , - K„x r . 


*nter line )f jet tube 


♦ , then: 


«« I h.-B./A, 


a’ ■-.•inkin’ ‘ in Vi a range ,f 0-jo0 degre -• i. 


■ ' i .i - ; 


tf. . i t i ? pshfu : to precisely 


o.--a . -,aV.. .v. . : . ; . vv;,h.. 














(3) Calculation of Eccentric Thrust 

Jet tubes, due to the effects of additional designated oscillations 
ami contained pressures, are capable of causing jet tube axis lines to 
give rise to radial displacement, creating thrust eccentricity. This is 
an important technical target of flexible connector heads. 

In initial configuration, the jet tube axis line equation is; 

X-X. _Y-Y._Z-Z. (-34) 

A. B. C. 


T n th * nations: 

X.tVfZt -are the coordinates of a certain point on the OZ 

a x l s . 

A,.B,,C. -is the directional number of the OZ axis, which is 

always, for X, = Y,=*o,A t — i?,=0,Z, and 0, equal to certain 

constants. If we let the thrust eccentricity be D< f and we solve 
el lit ions (35) and (34), we then have: 


A 


1 

1- 

Z-Z % 

A t B t C, 

v4, Bi C, 

V 

A . 

B, 

1 

+ 

5. C. 

5, C 

S 

+ 

C.A . 1 
C >4. 


n/ B\ + A] 


(35) 


(1) >\ngala~ Vslocity of Jet Tube Oscillation 

\r Ve t' 'i-> ^unction we already obtained for the angles of rotation 
arc.uni the various axes: 

• Ar)ia,-/(/C • bT)*6 = f{K • A T) /C = 0,1,.n 

A7“-equidist'.nt stapling intervals. 


— 0 ,-,)/AT 

a>y i = (a i -a,-,)/&T ( 36 ) 

G>. ( = (0,-0,-,)/Ar 


[»] 




The over-nil angular velocity is: 

( 37 ) 

(5) Stroke of activator tubes 

In the third section, we already solved for the instantaneous 
lengths of activator tubes L.t and *L,, with the zero position length 
of the activator* tube as 'L. . Fro a this, the stroke of the activator 
tubes is; 

r AL m , — L m , — L, ( 88 ) 

' ALyi^Lyi — L» 

(5) Spatial displacements caused by the effects of axial 
displacement compensation and contained jet tube pressures 

This compensation is automatically carried out when the command 
signal is zero. After that, one carries out measurements and 
calculations of its spatial displacement. This is also another important 
technical parameter of flexible connector heads- 

With jet tubes in their initial configuration ( pressure zero, 
command zero) the electrical zero value of the mechanical apparatus and 
the sensors should satisfy; 

After carrying out compensation for the effects of contained 
pressure, this should satisfy: 

\Xn— Ah;i<0 .l»|JV—We then have: 

Axial displacement: 

Radial displacement: 


a z = z,;-z,: } 

AX - XX > (89) 

AY-Y»-Y X \ ] 








Note: TJal t.s in all the equation i tni salon l.at i on s ire t ; mo*n In l >«: 

Length: The r ul L l r,,.r n a j i a. it, »n 15 a u-i; t,. All j t • 1 - ~ mi ts 

e.re H'i. 

Angles and angular velocities: these respectively use degree; and 
dig* ns par second as units. 

^ > r • * til ;)•* is Mrs: these respectively use kilogram and kilograms 
per square cm as units. 

Mo nent s of fore-*: V i 1 ogra ns/ net.er . 

Time: seconds. 
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